﻿#pragma once

#include <learn/alg/alg_projection.h>

namespace xi
{

namespace alg
{

struct MaxDistance
{
    double dist = 0;

    void accumulate(double d)
    {
        dist = std::max(dist, d);
    }

    double result() const
    {
        return dist;
    }
};

struct AverageDistance
{
    double dist = 0;
    std::size_t count = 0;

    void accumulate(double d)
    {
        dist += d;
        count++;
    }

    double result() const
    {
        return dist / count;
    }
};

struct MeanSquareDistance
{
    double dist = 0;
    std::size_t count = 0;

    void accumulate(double d)
    {
        dist += d * d;
        count++;
    }

    double result() const
    {
        return dist / count;
    }
};

template <std::size_t N> struct Alg_Distance
{
    using Vec = VecT<N, double>; // P

    /**
     * @brief 计算两个点集之间的 Hausdorff 距离
     *
     * @param[in] cloud1
     * @param[in] cloud2
     * @return double
     */
    static double distance(const std::vector<Vec> &cloud1, const std::vector<Vec> &cloud2)
    {
        double dist = 0;
        for (std::size_t i = 0; i < cloud1.size(); i++)
        {
            double d = std::numeric_limits<double>::max();
            for (std::size_t j = 0; j < cloud2.size(); j++)
                d = std::min(d, cloud1[i].distance(cloud2[j]));
            dist = std::max(dist, d);
        }
        return dist;
    }

    /**
     * @brief 计算采样点到曲线的 Hausdorff 距离
     *
     * @param[in] samples
     * @param[in] curve
     * @return double
     */
    template <class Metric = MaxDistance>
    static double distance(const std::vector<Vec> &samples, std::shared_ptr<Geom_Curve<N>> curve)
    {
        Metric m;
        for (std::size_t i = 0; i < samples.size(); i++)
        {
            double t = (curve->first_param() + curve->last_param()) / 2;
            t = Alg_Projection<N>::orthogonal2(curve, samples[i], t, Numeric::half);
            t = Alg_Projection<N>::nearest(curve, samples[i], t, Numeric::half);
            auto p = curve->value(t);
            m.accumulate(p.distance(samples[i]));
        }
        return m.result();
    }

    /**
     * @brief 计算两条曲线之间的 Hausdorff 距离
     *
     * @param[in] curve1
     * @param[in] curve2
     * @return double
     */
    static double distance(std::shared_ptr<Geom_Curve<N>> curve1, std::shared_ptr<Geom_Curve<N>> curve2)
    {
        static constexpr int n = 1000;

        auto step1 = (curve1->last_param() - curve1->first_param()) / (n - 1);
        auto step2 = (curve2->last_param() - curve2->first_param()) / (n - 1);
        std::vector<Vec> points1, points2;
        for (std::size_t i = 0; i < n; i++)
        {
            auto p = curve1->value(curve1->first_param() + i * step1);
            points1.push_back(p);
        }

        for (std::size_t i = 0; i < n; i++)
        {
            auto p = curve2->value(curve2->first_param() + i * step2);
            points2.push_back(p);
        }

        return distance(points1, points2);
    }

    /**
     * @brief 计算采样点到曲面的 Hausdorff 距离
     *
     * @param[in] samples
     * @param[in] surface
     * @return double
     */
    template <class Metric = MaxDistance>
    static double distance(const std::vector<Vec> &samples, std::shared_ptr<Geom_Surface<N>> surface)
    {
        Metric m;
        for (std::size_t i = 0; i < samples.size(); i++)
        {
            auto [u1, u2, v1, v2] = surface->bounds();
            double u = (u1 + u2) / 2, v = (v1 + v2) / 2;
            auto uv = Alg_Projection<N>::torus(surface, samples[i], u, v, Numeric::half);
            auto p = surface->value(uv.x(), uv.y());
            m.accumulate(p.distance(samples[i]));
        }
        return m.result();
    }

    /**
     * @brief 计算采样点到曲面的 Hausdorff 距离
     *
     * @param[in] samples
     * @param[in] params
     * @param[in] surface
     * @return double
     */
    template <class Metric = MaxDistance>
    static double distance(const std::vector<Vec> &samples, const std::vector<Vec2> &params,
                           std::shared_ptr<Geom_Surface<N>> surface)
    {
        Metric m;
        for (std::size_t i = 0; i < samples.size(); i++)
        {
            double u = params[i].x(), v = params[i].y();
            auto p = surface->value(u, v);
            m.accumulate(p.distance(samples[i]));
        }
        return m.result();
    }
};

} // namespace alg

} // namespace xi
